*
* $Id$
*
* $Log: gheish.F,v $
* Revision 1.1.1.1  2002/07/24 15:56:27  rdm
* initial import into CVS
*
* Revision 1.1.1.1  2002/06/16 15:18:40  hristov
* Separate distribution  of Geant3
*
* Revision 1.1.1.1  1999/05/18 15:55:19  fca
* AliRoot sources
*
* Revision 1.1.1.1  1995/10/24 10:21:15  cernlib
* Geant
*
*
#include "geant321/pilot.h"
*CMZ :  3.21/02 29/03/94  15.41.38  by  S.Giani
*-- Author :
      SUBROUTINE GHEISH
C
C *** MAIN STEERING FOR HADRON SHOWER DEVELOPMENT ***
C *** NVE 15-JUN-1988 CERN GENEVA ***
C
C CALLED BY : GUHADR (USER ROUTINE)
C ORIGIN : F.CARMINATI, H.FESEFELDT
C                       ROUTINES : CALIM  16-SEP-1987
C                                  SETRES 19-AUG-1985
C                                  INTACT 06-OCT-1987
C
#include "geant321/gcbank.inc"
#include "geant321/gcjloc.inc"
#include "geant321/gccuts.inc"
#include "geant321/gcflag.inc"
#include "geant321/gckine.inc"
#include "geant321/gcking.inc"
#include "geant321/gcmate.inc"
#include "geant321/gcphys.inc"
#include "geant321/gctmed.inc"
#include "geant321/gctrak.inc"
#include "geant321/gsecti.inc"
#include "geant321/gcunit.inc"
C --- GHEISHA COMMONS ---
#include "geant321/mxgkgh.inc"
#include "geant321/s_blankp.inc"
#include "geant321/s_consts.inc"
#include "geant321/s_event.inc"
#include "geant321/s_prntfl.inc"
C
C --- "NEVENT" CHANGED TO "KEVENT" IN COMMON /CURPAR/ DUE TO CLASH ---
C --- WITH VARIABLE "NEVENT" IN GEANT COMMON ---
C
      PARAMETER (MXGKCU=MXGKGH)
      COMMON /CURPAR /WEIGHT(10),DDELTN,IFILE,IRUN,NEVT,KEVENT,SHFLAG,
     $                ITHST,ITTOT,ITLST,IFRND,TOFCUT,CMOM(5),CENG(5),
     $                RS,S,ENP(10),NP,NM,NN,NR,NO,NZ,IPA(MXGKCU),
     $                ATNO2,ZNO2
C
C --- "IPART" CHANGED TO "KPART" IN COMMON /RESULT/ DUE TO CLASH ---
C --- WITH VARIABLE "IPART" IN GEANT COMMON ---
C
      COMMON /RESULT/ XEND,YEND,ZEND,RCA,RCE,AMAS,NCH,TOF,PX,PY,PZ,
     $                USERW,INTCT,P,EN,EK,AMASQ,DELTN,ITK,NTK,KPART,IND,
     $                LCALO,ICEL,SINL,COSL,SINP,COSP,
     $                XOLD,YOLD,ZOLD,POLD,PXOLD,PYOLD,PZOLD,
     $                XSCAT,YSCAT,ZSCAT,PSCAT,PXSCAT,PYSCAT,PZSCAT
                      REAL NCH,INTCT
C
C --- "ABSL(21)" CHANGED TO "ABSLTH(21)" IN COMMON /MAT/ DUE TO CLASH ---
C --- WITH VARIABLE "ABSL" IN GEANT COMMON ---
C
      COMMON /MAT/ LMAT,
     $             DEN(21),RADLTH(21),ATNO(21),ZNO(21),ABSLTH(21),
     $             CDEN(21),MDEN(21),X0DEN(21),X1DEN(21),RION(21),
     $             MATID(21),MATID1(21,24),PARMAT(21,10),
     $             IFRAT,IFRAC(21),FRAC1(21,10),DEN1(21,10),
     $             ATNO1(21,10),ZNO1(21,10)
C
      DIMENSION IPELOS(40)
      SAVE IDEOL
C
C --- TRANSFER GEANT CUT-OFFS INTO GHEISHA VALUES ---
      DIMENSION CUTS(5)
      EQUIVALENCE (CUTS(1),CUTGAM)
      DIMENSION RNDM(1)
C
#include "geant321/pcodim.inc"
#include "geant321/pcodat.inc"
C
C --- DENOTE STABLE PARTICLES ACCORDING TO GHEISHA CODE ---
C --- STABLE : GAMMA, NEUTRINO, ELECTRON, PROTON AND HEAVY FRAGMENTS ---
C --- WHEN STOPPING THESE PARTICLES ONLY LOOSE THEIR KINETIC ENERGY ---
      DATA IPELOS/
     $             1,   1,   0,   1,   0,   0,   0,   0,
     $             0,   0,   0,   0,   0,   1,   0,   0,
     $             0,   0,   0,   0,   0,   0,   0,   0,
     $             0,   0,   0,   0,   0,   1,   1,   1,
     $             0,   0,   1,   1,   1,   1,   1,   1/
C
C --- LOWERBOUND OF KINETIC ENERGY BIN IN N CROSS-SECTION TABLES ---
      DATA TEKLOW /0.0001/
C
C --- KINETIC ENERGY TO SWITCH FROM "CASN" TO "GNSLWD" FOR N CASCADE ---
      DATA SWTEKN /0.05/
C
      DATA IDEOL/0/
C
C --- INITIALIZE RELEVANT GHEISHA VARIABLES IN CASE NOT DONE ALREADY ---
      IF (IFINIT(4) .EQ. 0) CALL GHEINI
C
C --- SET THE INTERACTION MECHANISM TO "HADR" ---
      KCASE=NAMEC(12)
C
C --- SET GHEISHA PRINTING FLAGS ACCORDING TO "DEBUG" STEERING CARD --
      IF (IDEOL .EQ. IDEBUG) GO TO 9000
C
      IF (IDEBUG .NE. 1) GO TO 9001
C
C --- SET SELECTED DEBUGGING FLAGS ---
      DO 9002 LL=1,10
      IF ((ISWIT(LL) .LE. 100) .OR. (ISWIT(LL) .GT. 110)) GO TO 9002
      JJ=ISWIT(LL)-100
      NPRT(JJ)=.TRUE.
 9002 CONTINUE
      GO TO 9000
C
C --- NO DEBUGGING SELECTED ---
 9001 CONTINUE
      DO 9003 LL=1,10
      NPRT(LL)=.FALSE.
 9003 CONTINUE
      IDEOL=IDEBUG
C
 9000 CONTINUE
C
C --- SET THE GHEISHA PARTICLE TYPE TO THE ONE OF GEANT ---
      IF(IPART.GT.79) THEN
         IF(ISTOP.EQ.0) GOTO 9999
         JPA = LQ(JPART-IPART)
         AMAS=Q(JPA+7)
         NCH =Q(JPA+8)
         KPART=-IPART
         GOTO 107
      ENDIF
      NETEST=IKPART(KPART)
      IF ((NETEST .EQ. IPART) .OR. (ISTOP .NE. 0)) GO TO 9004
C
      PRINT 8881,IPART,KPART,ISTOP
 8881 FORMAT(' *GHEISH* IPART,KPART = ',2(I3,1X),' ISTOP = ',I3/
     $ ' *GHEISH* ======> PARTICLE TYPES DO NOT MATCH <=======')
      STOP
C
 9004 CONTINUE
      KPART=KIPART(IPART)
      KKPART=KPART
      AMAS=RMASS(KPART)
      NCH=RCHARG(KPART)
C
C --- TRANSPORT THE TRACK NUMBER TO GHEISHA AND INITIALISE SOME NUMBERS
 107  NTK=ITRA
      INTCT=0.0
      NEXT=1
      NTOT=0
      TOF=0.0
C
C --- FILL RESULT COMMON FOR THIS TRACK WITH GEANT VALUES ---
C --- CALIM CODE ---
      XEND=VECT(1)
      YEND=VECT(2)
      ZEND=VECT(3)
      PX=VECT(4)
      PY=VECT(5)
      PZ=VECT(6)
      USERW=UPWGHT
C --- SETRES CODE ---
      P=VECT(7)
      AMASQ=AMAS*AMAS
      EN=SQRT(AMASQ+P*P)
      EK=ABS(EN-ABS(AMAS))
      ENOLD=EN
C
      SINL=0.0
      COSL=1.0
      SINP=0.0
      COSP=1.0
C
      IF (ABS(P) .LE. 1.0E-10) GO TO 1
      SINL=PZ
      COSL=SQRT(ABS(1.0-SINL**2))
C
 1    CONTINUE
      CALL GRNDM(RNDM,1)
      PHI=RNDM(1)*TWPI
      IF ((PX .EQ. 0.0) .AND. (PY .EQ. 0.0)) GOTO 3
      IF (ABS(PX) .LT. 1.E-10) GOTO 2
      PHI=ATAN2(PY,PX)
      GOTO 3
C
 2    CONTINUE
      IF (PY .GT. 0.0) PHI=PI/2.0
      IF (PY .LE. 0.0) PHI=3.0*PI/2.0
C
 3    CONTINUE
      SINP=SIN(PHI)
      COSP=COS(PHI)
C
C --- SET GHEISHA INDEX FOR THE CURRENT MEDIUM ALWAYS TO 1 ---
      IND=1
C
C --- TRANSFER GLOBAL MATERIAL CONSTANTS FOR CURRENT MEDIUM ---
C --- DETAILED DATA FOR COMPOUNDS IS OBTAINED VIA ROUTINE COMPO ---
      ATNO(IND+1)=A
      ZNO(IND+1)=Z
      DEN(IND+1)=DENS
      RADLTH(IND+1)=RADL
      ABSLTH(IND+1)=ABSL
C
C --- SETUP PARMAT FOR PHYSICS STEERING ---
      PARMAT(IND+1,5)=0.0
      PARMAT(IND+1,8)=IPFIS
      PARMAT(IND+1,9)=0.0
      PARMAT(IND+1,10)=0.0
      JTMN=LQ(JTM)
      IF (JTMN .LE. 0) GO TO 4
      PARMAT(IND+1,5)=Q(JTMN+26)
 4    CONTINUE
C
C --- CHECK WHETHER PARTICLE IS STOPPING OR NOT ---
      IF (ISTOP .EQ. 0) GO TO 5
C
      IF (NPRT(9)) PRINT 1000,KPART
 1000 FORMAT(' *GHEISH* STOPPING GHEISHA PARTICLE ',I3)
      CALL GHSTOP
C --- IN CASE OF DECAY OF PARTICLE OR USER PARTICLE ==> RETURN ---
      IF (LMEC(NMEC) .EQ. 5 .OR. KPART .LT. 0) GO TO 9999
C --- IN CASE OF HAD. INT. WITH GENERATION OF SEC. ==> GO TO 40 ---
      IF (IHADR .NE. 2) GO TO 40
C --- ALSO DEPOSIT REST MASS ENERGY FOR IN-STABLE PARTICLES ---
      IF (IPELOS(KPART) .EQ. 0) DESTEP=DESTEP+ABS(RMASS(KPART))
      GO TO 9999
  5   CONTINUE
C
C --- INDICATE LIGHT (<= PI) AND HEAVY PARTICLES (HISTORICALLY) ---
C --- CALIM CODE ---
      J=2
      TEST=RMASS(7)-0.001
      IF (ABS(AMAS) .LT. TEST) J=1
C
C *** DIVISION INTO VARIOUS INTERACTION CHANNELS DENOTED BY "INT" ***
C THE CONVENTION FOR "INT" IS THE FOLLOWING
C
C INT  = -1 REACTION CROSS SECTIONS NOT YET TABULATED/PROGRAMMED
C      =  0 NO INTERACTION
C      =  1 ELEASTIC SCATTERING
C      =  2 INELASTIC SCATTERING
C      =  3 NUCLEAR FISSION WITH INELEASTIC SCATTERING
C      =  4 NEUTRON CAPTURE
C
C --- INTACT CODE ---
      KK=ABS(Q(JMA+11))
      ALAM1=0.0
      CALL GRNDM(RNDM,1)
      RAT=RNDM(1)*ALAM
      NMEC=NMEC+1
      ATNO2=A
      ZNO2 =Z
C
      DO 6 K=1,KK
      IF (KK .LE. 0) GO TO 6
C
      IF (KK .EQ. 1) GO TO 7
      ATNO2=Q(JMIXT+K)
      ZNO2 =Q(JMIXT+K+KK)
C
 7    CONTINUE
C
C --- TRY FOR ELASTIC SCATTERING ---
      INT=1
      LMEC(NMEC)=13
      ALAM1=ALAM1+AIEL(K)
      IF (RAT .LT. ALAM1) GO TO 8
C
C --- TRY FOR INELASTIC SCATTERING ---
      INT=2
      LMEC(NMEC)=20
      ALAM1=ALAM1+AIIN(K)
      IF (RAT .LT. ALAM1) GO TO 8
C
C --- TRY FOR NUCLEAR FISSION WITH INELASTIC SCATTERING ---
      INT=3
      LMEC(NMEC)=15
      ALAM1=ALAM1+AIFI(K)
      IF (RAT .LT. ALAM1) GO TO 8
C
C --- TRY FOR NEUTRON CAPTURE ---
      INT=4
      LMEC(NMEC)=18
      ALAM1=ALAM1+AICA(K)
      IF (RAT .LT. ALAM1) GO TO 8
C
 6    CONTINUE
C --- NO REACTION SELECTED ==> ELASTIC SCATTERING ---
      INT=1
      LMEC(NMEC)=13
C
C *** TAKE ACTION ACCORDING TO SELECTED REACTION CHANNEL ***
C --- FOLLOWING CODE IS A TRANSLATION OF "CALIM" INTO GEANT JARGON ---
C
 8    CONTINUE
      IF (NPRT(9)) PRINT 1001,INT
 1001 FORMAT(' *GHEISH* INTERACTION TYPE CHOSEN INT = ',I3)
C
C --- IN CASE OF NO INTERACTION OR UNKNOWN CROSS SECTIONS ==> DONE ---
      IF (INT .LE. 0) GO TO 40
C
C --- IN CASE OF NON-ELASTIC SCATTERING AND NO GENERATION OF SEC. ---
C --- PARTICLES DEPOSIT TOTAL PARTICLE ENERGY AND RETURN ---
      IF ((INT .EQ. 1) .OR. (IHADR .NE. 2)) GO TO 9
      ISTOP=2
      DESTEP=DESTEP+EN
      NGKINE=0
      GO TO 9999
C
 9    CONTINUE
      IF (INT .NE. 4) GO TO 10
C
C --- NEUTRON CAPTURE ---
      IF (NPRT(9)) PRINT 2000
 2000 FORMAT(' *GHEISH* ROUTINE CAPTUR WILL BE CALLED')
      ISTOP=1
      CALL CAPTUR(NOPT)
      GO TO 40
C
 10   CONTINUE
      IF (INT .NE. 3) GO TO 11
C --- NUCLEAR FISSION ---
      IF (NPRT(9)) PRINT 2001
 2001 FORMAT(' *GHEISH* ROUTINE FISSIO WILL BE CALLED')
      ISTOP=1
      TKIN=FISSIO(EK)
      GO TO 40
C
 11   CONTINUE
C
C --- ELASTIC AND INELASTIC SCATTERING ---
      PV( 1,MXGKPV)=P*PX
      PV( 2,MXGKPV)=P*PY
      PV( 3,MXGKPV)=P*PZ
      PV( 4,MXGKPV)=EN
      PV( 5,MXGKPV)=AMAS
      PV( 6,MXGKPV)=NCH
      PV( 7,MXGKPV)=TOF
      PV( 8,MXGKPV)=KPART
      PV( 9,MXGKPV)=0.
      PV(10,MXGKPV)=USERW
C
C --- ADDITIONAL PARAMETERS TO SIMULATE FERMI MOTION AND EVAPORATION ---
      DO 111 JENP=1,10
         ENP(JENP)=0.
 111  CONTINUE
      ENP(5)=EK
      ENP(6)=EN
      ENP(7)=P
C
      IF (INT .NE. 1) GO TO 12
C
C *** ELASTIC SCATTERING PROCESSES ***
C
C --- ONLY NUCLEAR INTERACTIONS FOR HEAVY FRAGMENTS ---
C      IF ((KPART .GE. 30) .AND. (KPART .LE. 32)) GO TO 35
      IF (((KPART.GE.30).AND.(KPART.LE.32)).OR.
     $    ((KPART.GE.36).AND.(KPART.LE.40))) GO TO 35

C
C --- NORMAL ELASTIC SCATTERING FOR LIGHT MEDIA ---
      IF (ATNO2 .LT. 1.5) GO TO 35
C
C --- COHERENT ELASTIC SCATTERING FOR HEAVY MEDIA ---
      IF (NPRT(9)) PRINT 2002
 2002 FORMAT(' *GHEISH* ROUTINE COSCAT WILL BE CALLED')
      CALL COSCAT
      GO TO 40
C
C *** NON-ELASTIC SCATTERING PROCESSES ***
 12   CONTINUE
C
C --- ONLY NUCLEAR INTERACTIONS FOR HEAVY FRAGMENTS ---
C      IF ((KPART .GE. 30) .AND. (KPART .LE. 32)) GO TO 35
      IF (((KPART.GE.30).AND.(KPART.LE.32)).OR.
     $    ((KPART.GE.36).AND.(KPART.LE.40))) GO TO 35

C
C *** USE SOMETIMES NUCLEAR REACTION ROUTINE "NUCREC" FOR LOW ENERGY ***
C *** PROTON AND NEUTRON SCATTERING ***
      CALL GRNDM(RNDM,1)
      TEST1=RNDM(1)
      TEST2=4.5*(EK-0.01)
      IF ((KPART .EQ. 14) .AND. (TEST1 .GT. TEST2)) GO TO 85
      IF ((KPART .EQ. 16) .AND. (TEST1 .GT. TEST2)) GO TO 86
C
C *** FERMI MOTION AND EVAPORATION ***
      TKIN=CINEMA(EK)
      PV( 9,MXGKPV)=TKIN
      ENP(5)=EK+TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=(ENP(6)-AMAS)*(ENP(6)+AMAS)
      ENP(7)=SQRT(ABS(ENP(7)))
      TKIN=FERMI(ENP(5))
      ENP(5)=ENP(5)+TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=(ENP(6)-AMAS)*(ENP(6)+AMAS)
      ENP(7)=SQRT(ABS(ENP(7)))
      TKIN=EXNU(ENP(5))
      ENP(5)=ENP(5)-TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=(ENP(6)-AMAS)*(ENP(6)+AMAS)
      ENP(7)=SQRT(ABS(ENP(7)))
C
C *** IN CASE OF ENERGY ABOVE CUT-OFF LET THE PARTICLE CASCADE ***
      TEST=ABS(CHARGE)
      IF ((TEST .GT. 1.0E-10) .AND. (ENP(5) .GT. CUTHAD)) GO TO 35
      IF ((TEST .LE. 1.0E-10) .AND. (ENP(5) .GT. CUTNEU)) GO TO 35
C
C --- SECOND CHANCE FOR ANTI-BARYONS DUE TO POSSIBLE ANNIHILATION ---
      IF ((AMAS .GE. 0.0) .OR. (KPART .LE. 14)) GO TO 13
      ANNI=1.3*P
      IF (ANNI .GT. 0.4) ANNI=0.4
      CALL GRNDM(RNDM,1)
      TEST=RNDM(1)
      IF (TEST .GT. ANNI) GO TO 35
C
C *** PARTICLE WITH ENERGY BELOW CUT-OFF ***
C --- ==> ONLY NUCLEAR EVAPORATION AND QUASI-ELASTIC SCATTERING ---
 13   CONTINUE
C
      ISTOP=3
C
      IF (NPRT(9)) PRINT 1002,KPART,EK,EN,P,ENP(5),ENP(6),ENP(7)
 1002 FORMAT(' *GHEISH* ENERGY BELOW CUT-OFF FOR GHEISHA PARTICLE ',I3/
     $ ' EK,EN,P,ENP(5),ENP(6),ENP(7) = ',6(G12.5,1X))
C
      IF ((KPART .NE. 14) .AND. (KPART .NE. 16)) GO TO 14
      IF (KPART .EQ. 16) GO TO 86
C
C --- SLOW PROTON ---
 85   CONTINUE
      IF (NPRT(9)) PRINT 2003,EK,KPART
 2003 FORMAT(' *GHEISH* ROUTINE NUCREC WILL BE CALLED',
     $ ' EK = ',G12.5,' GEV  KPART = ',I3)
      CALL NUCREC(NOPT,2)
C
      IF (NOPT .NE. 0) GO TO 50
C
      IF (NPRT(9)) PRINT 2004,EK,KPART
 2004 FORMAT(' *GHEISH* ROUTINE COSCAT WILL BE CALLED',
     $ ' EK = ',G12.5,' GEV  KPART = ',I3)
      CALL COSCAT
      GO TO 40
C
C --- SLOW NEUTRON ---
 86   CONTINUE
      IF (NPRT(9)) PRINT 2015
      NUCFLG=0
      CALL GNSLWD(NUCFLG,INT,NFL,TEKLOW)
      IF (NUCFLG .NE. 0) GO TO 50
      GO TO 40
C
C --- OTHER SLOW PARTICLES ---
 14   CONTINUE
      IPA(1)=KPART
C --- DECIDE FOR PROTON OR NEUTRON TARGET ---
      IPA(2)=16
      CALL GRNDM(RNDM,1)
      TEST1=RNDM(1)
      TEST2=ZNO2/ATNO2
      IF (TEST1 .LT. TEST2) IPA(2)=14
      AVERN=0.0
      NFL=1
      IF (IPA(2) .EQ. 16) NFL=2
      IPPP=KPART
      IF (NPRT(9)) PRINT 2005
 2005 FORMAT(' *GHEISH* ROUTINE TWOB WILL BE CALLED')
      CALL TWOB(IPPP,NFL,AVERN)
      GOTO 40
C
C --- INITIALISATION OF CASCADE QUANTITIES ---
 35   CONTINUE
C
C *** CASCADE GENERATION ***
C --- CALCULATE FINAL STATE MULTIPLICITY AND LONGITUDINAL AND ---
C --- TRANSVERSE MOMENTUM DISTRIBUTIONS ---
C
C --- FIXED PARTICLE TYPE TO STEER THE CASCADE ---
      KKPART=KPART
C
C --- NO CASCADE FOR LEPTONS ---
      IF (KKPART .LE. 6) GO TO 9999
C
C *** WHAT TO DO WITH "NEW PARTICLES" FOR GHEISHA ?????? ***
C --- RETURN FOR THE TIME BEING ---
C      IF (KKPART .GE. 35) GO TO 9999
      IF (KKPART .EQ. 35) GO TO 9999
C
C --- CASCADE OF HEAVY FRAGMENTS
C      IF ((KKPART .GE. 30) .AND. (KKPART .LE. 32)) GO TO 390
      IF (((KKPART.GE.30).AND.(KKPART.LE.32)).OR.
     $    ((KKPART.GE.36).AND.(KKPART.LE.40))) GO TO 390

C
C --- INITIALIZE THE IPA ARRAY ---
      CALL VZERO(IPA(1),MXGKCU)
C
C --- CASCADE OF OMEGA - AND OMEGA - BAR ---
      IF (KKPART .EQ. 33) GO TO 330
      IF (KKPART .EQ. 34) GO TO 331
C
      NVEPAR=KKPART-17
      IF (NVEPAR .LE. 0) GO TO 15
      GO TO (318,319,320,321,322,323,324,325,326,327,328,329),NVEPAR
C
 15   CONTINUE
      NVEPAR=KKPART-6
      GO TO (307,308,309,310,311,312,313,314,315,316,317,318),NVEPAR
C
C --- PI+ CASCADE ---
 307  CONTINUE
      IF (NPRT(9)) PRINT 2006
 2006 FORMAT(' *GHEISH* ROUTINE CASPIP WILL BE CALLED')
      CALL CASPIP(J,INT,NFL)
      GO TO 40
C
C --- PI0 ==> NO CASCADE ---
 308  CONTINUE
      GO TO 40
C
C --- PI- CASCADE ---
 309  CONTINUE
      IF (NPRT(9)) PRINT 2007
 2007 FORMAT(' *GHEISH* ROUTINE CASPIM WILL BE CALLED')
      CALL CASPIM(J,INT,NFL)
      GO TO 40
C
C --- K+ CASCADE ---
 310  CONTINUE
      IF (NPRT(9)) PRINT 2008
 2008 FORMAT(' *GHEISH* ROUTINE CASKP WILL BE CALLED')
      CALL CASKP(J,INT,NFL)
      GO TO 40
C
C --- K0 CASCADE ---
 311  CONTINUE
      IF (NPRT(9)) PRINT 2009
 2009 FORMAT(' *GHEISH* ROUTINE CASK0 WILL BE CALLED')
      CALL CASK0(J,INT,NFL)
      GO TO 40
C
C --- K0 BAR CASCADE ---
 312  CONTINUE
      IF (NPRT(9)) PRINT 2010
 2010 FORMAT(' *GHEISH* ROUTINE CASK0B WILL BE CALLED')
      CALL CASK0B(J,INT,NFL)
      GO TO 40
C
C --- K- CASCADE ---
 313  CONTINUE
      IF (NPRT(9)) PRINT 2011
 2011 FORMAT(' *GHEISH* ROUTINE CASKM WILL BE CALLED')
      CALL CASKM(J,INT,NFL)
      GO TO 40
C
C --- PROTON CASCADE ---
 314  CONTINUE
      IF (NPRT(9)) PRINT 2012
 2012 FORMAT(' *GHEISH* ROUTINE CASP WILL BE CALLED')
      CALL CASP(J,INT,NFL)
      GO TO 40
C
C --- PROTON BAR CASCADE ---
 315  CONTINUE
      IF (NPRT(9)) PRINT 2013
 2013 FORMAT(' *GHEISH* ROUTINE CASPB WILL BE CALLED')
      CALL CASPB(J,INT,NFL)
      GO TO 40
C
C --- NEUTRON CASCADE ---
 316  CONTINUE
      NUCFLG=0
      IF (EK .GT. SWTEKN) THEN
         CALL CASN(J,INT,NFL)
         IF (NPRT(9)) PRINT 2014
 2014 FORMAT(' *GHEISH* ROUTINE CASN WILL BE CALLED')
      ELSE
         CALL GNSLWD(NUCFLG,INT,NFL,TEKLOW)
         IF (NPRT(9)) PRINT 2015
 2015 FORMAT(' *GHEISH* ROUTINE GNSLWD WILL BE CALLED')
      ENDIF
      IF (NUCFLG .NE. 0) GO TO 50
      GO TO 40
C
C --- NEUTRON BAR CASCADE ---
 317  CONTINUE
      IF (NPRT(9)) PRINT 2016
 2016 FORMAT(' *GHEISH* ROUTINE CASNB WILL BE CALLED')
      CALL CASNB(J,INT,NFL)
      GO TO 40
C
C --- LAMBDA CASCADE ---
 318  CONTINUE
      IF (NPRT(9)) PRINT 2017
 2017 FORMAT(' *GHEISH* ROUTINE CASL0 WILL BE CALLED')
      CALL CASL0(J,INT,NFL)
      GO TO 40
C
C --- LAMBDA BAR CASCADE ---
 319  CONTINUE
      IF (NPRT(9)) PRINT 2018
 2018 FORMAT(' *GHEISH* ROUTINE CASAL0 WILL BE CALLED')
      CALL CASAL0(J,INT,NFL)
      GO TO 40
C
C --- SIGMA + CASCADE ---
 320  CONTINUE
      IF (NPRT(9)) PRINT 2019
 2019 FORMAT(' *GHEISH* ROUTINE CASSP WILL BE CALLED')
      CALL CASSP(J,INT,NFL)
      GO TO 40
C
C --- SIGMA 0 ==> NO CASCADE ---
 321  CONTINUE
      GO TO 40
C
C --- SIGMA - CASCADE ---
 322  CONTINUE
      IF (NPRT(9)) PRINT 2020
 2020 FORMAT(' *GHEISH* ROUTINE CASSM WILL BE CALLED')
      CALL CASSM(J,INT,NFL)
      GO TO 40
C
C --- SIGMA + BAR CASCADE ---
 323  CONTINUE
      IF (NPRT(9)) PRINT 2021
 2021 FORMAT(' *GHEISH* ROUTINE CASASP WILL BE CALLED')
      CALL CASASP(J,INT,NFL)
      GO TO 40
C
C --- SIGMA 0 BAR ==> NO CASCADE ---
 324  CONTINUE
      GO TO 40
C
C --- SIGMA - BAR CASCADE ---
 325  CONTINUE
      IF (NPRT(9)) PRINT 2022
 2022 FORMAT(' *GHEISH* ROUTINE CASASM WILL BE CALLED')
      CALL CASASM(J,INT,NFL)
      GO TO 40
C
C --- XI 0 CASCADE ---
 326  CONTINUE
      IF (NPRT(9)) PRINT 2023
 2023 FORMAT(' *GHEISH* ROUTINE CASX0 WILL BE CALLED')
      CALL CASX0(J,INT,NFL)
      GO TO 40
C
C --- XI - CASCADE ---
 327  CONTINUE
      IF (NPRT(9)) PRINT 2024
 2024 FORMAT(' *GHEISH* ROUTINE CASXM WILL BE CALLED')
      CALL CASXM(J,INT,NFL)
      GO TO 40
C
C --- XI 0 BAR CASCADE ---
 328  CONTINUE
      IF (NPRT(9)) PRINT 2025
 2025 FORMAT(' *GHEISH* ROUTINE CASAX0 WILL BE CALLED')
      CALL CASAX0(J,INT,NFL)
      GO TO 40
C
C --- XI - BAR CASCADE ---
 329  CONTINUE
      IF (NPRT(9)) PRINT 2026
 2026 FORMAT(' *GHEISH* ROUTINE CASAXM WILL BE CALLED')
      CALL CASAXM(J,INT,NFL)
      GO TO 40
C
C --- OMEGA - CASCADE ---
 330  CONTINUE
      IF (NPRT(9)) PRINT 2027
 2027 FORMAT(' *GHEISH* ROUTINE CASOM WILL BE CALLED')
      CALL CASOM(J,INT,NFL)
      GO TO 40
C
C --- OMEGA - BAR CASCADE ---
 331  CONTINUE
      IF (NPRT(9)) PRINT 2028
 2028 FORMAT(' *GHEISH* ROUTINE CASAOM WILL BE CALLED')
      CALL CASAOM(J,INT,NFL)
      GO TO 40
C
C --- HEAVY FRAGMENT CASCADE ---
 390  CONTINUE
      IF (NPRT(9)) PRINT 2090
 2090 FORMAT(' *GHEISH* ROUTINE CASFRG WILL BE CALLED')
      NUCFLG=0
      CALL CASFRG(NUCFLG,INT,NFL)
      IF (NUCFLG .NE. 0) GO TO 50
C
C *** CHECK WHETHER THERE ARE NEW PARTICLES GENERATED ***
 40   CONTINUE
      IF ((NTOT .NE. 0) .OR. (KKPART .NE. KPART)) GO TO 50
C
C --- NO SECONDARIES GENERATED AND PARTICLE IS STILL THE SAME ---
C --- ==> COPY EVERYTHING BACK IN THE CURRENT GEANT STACK ---
      NGKINE=0
      TOFG=TOFG+TOF*0.5E-10
C --- In case of crazy momentum value ==> no change to GEANT stack ---
      IF (P .LT. 0.) GO TO 41
      VECT(4)=PX
      VECT(5)=PY
      VECT(6)=PZ
      VECT(7)=P
      GETOT=EN
      GEKIN=EK
C --- CHECK KINETIC ENERGY ---
      CALL G3EKBIN
      EDEP=ABS(ENOLD-EN)
      RMASSI=EN-EK
      IF (NPRT(9) .AND. (EN .GT. ENOLD))
     $ PRINT 8888,EDEP,ENOLD,EN,EK,RMASSI
 8888 FORMAT(' *GHEISH* EDEP,ENOLD,EN,EK,M = ',5(G12.5,1X)/
     $ ' *GHEISH* =======> EDEP WOULD BE NEGATIVE <========')
      IF (ISTOP .EQ. 0) DESTEP=DESTEP+EDEP
C
C --- RE-INITIALIZE THE PROBABILITY FOR HADRONIC INTERACTION ---
 41   CONTINUE
      CALL GRNDM(RNDM,1)
      IF ((RNDM(1) .LE. 0.) .OR. (RNDM(1) .GE. 1.)) GO TO 41
      ZINTHA=-LOG(RNDM(1))
      SLHADR=SLENG
      STEPHA=1.0E10
C
      NVEDUM=KIPART(IPART)
      IF (NPRT(9)) PRINT 1003,NTOT,IPART,KPART,KKPART,NVEDUM
 1003 FORMAT(' *GHEISH* NO SEC. GEN. NTOT,IPART,KPART,KKPART,KIPART = ',
     $ 5(I3,1X)/
     $ ' CURRENT PARTICLE ON THE STACK AGAIN')
      GO TO 9999
C
C *** CURRENT PARTICLE IS NOT THE SAME AS IN THE BEGINNING OR/AND ***
C *** ONE OR MORE SECONDARIES HAVE BEEN GENERATED ***
 50   CONTINUE
C
      NVEDUM=KIPART(IPART)
      IF (NPRT(9)) PRINT 1004,NTOT,IPART,KPART,KKPART,NVEDUM
 1004 FORMAT(' *GHEISH* SEC. GEN. NTOT,IPART,KPART,KKPART,KIPART = ',
     $ 5(I3,1X))
C
C --- INITIAL PARTICLE TYPE HAS BEEN CHANGED ==> PUT NEW TYPE ON ---
C --- THE GEANT TEMPORARY STACK ---
C
C --- MAKE CHOICE BETWEEN K0 LONG / K0 SHORT ---
      IF ((KPART .NE. 11) .AND. (KPART .NE. 12)) GO TO 52
      CALL GRNDM(RNDM,1)
      KPART=11.5+RNDM(1)
C
 52   CONTINUE
      ITY=IKPART(KPART)
      LNVE=LQ(JPART-ITY)
      IF (LNVE .LE. 0) PRINT 1234,NTOT,ITY,LNVE
 1234 FORMAT('0*GHEISH* 1234 NTOT,ITY,LNVE = ',3(I10,1X))
      IF (LNVE .LE. 0) STOP
      IF (ISTOP .EQ. 0) ISTOP=1
C
C --- IN CASE THE NEW PARTICLE IS A NEUTRINO ==> FORGET IT ---
      IF (KPART .EQ. 2) GO TO 60
C
C --- PUT PARTICLE ON THE STACK ---
      GKIN(1,1)=PX*P
      GKIN(2,1)=PY*P
      GKIN(3,1)=PZ*P
      GKIN(4,1)=SQRT(P*P+RMASS(KPART)**2)
      GKIN(5,1)=ITY
      TOFD(1)=TOF*0.5E-10
      NGKINE = 1
      GPOS(1,1) = VECT(1)
      GPOS(2,1) = VECT(2)
      GPOS(3,1) = VECT(3)
C
      IF (NPRT(9)) PRINT 1005,ITY,NGKINE
 1005 FORMAT(' *GHEISH* GEANT PART. ',I3,' PUT ONTO STACK AT POS. ',I3)
C
C *** CHECK WHETHER SECONDARIES HAVE BEEN GENERATED AND COPY THEM ***
C *** ALSO ON THE GEANT STACK ***
 60   CONTINUE
C
C --- ALL QUANTITIES ARE TAKEN FROM THE GHEISHA STACK WHERE THE ---
C --- CONVENTION IS THE FOLLOWING ---
C
C EVE(INDEX+ 1)= X
C EVE(INDEX+ 2)= Y
C EVE(INDEX+ 3)= Z
C EVE(INDEX+ 4)= NCAL
C EVE(INDEX+ 5)= NCELL
C EVE(INDEX+ 6)= MASS
C EVE(INDEX+ 7)= CHARGE
C EVE(INDEX+ 8)= TOF
C EVE(INDEX+ 9)= PX
C EVE(INDEX+10)= PY
C EVE(INDEX+11)= PZ
C EVE(INDEX+12)= TYPE
C
      IF (NTOT .LE. 0) GO TO 9999
C
C --- ONE OR MORE SECONDARIES HAVE BEEN GENERATED ---
      DO 61 L=1,NTOT
      INDEX=(L-1)*12
      JND=EVE(INDEX+12)
C
C --- MAKE CHOICE BETWEEN K0 LONG / K0 SHORT ---
      IF ((JND .NE. 11) .AND. (JND .NE. 12)) GO TO 63
      CALL GRNDM(RNDM,1)
      JND=11.5+RNDM(1)
C
C --- FORGET ABOUT NEUTRINOS ---
 63   CONTINUE
      IF (JND .EQ. 2) GO TO 61
C
C --- SWITH TO GEANT QUANTITIES ---
      ITY=IKPART(JND)
      JTY=LQ(JPART-ITY)
      IF (JTY .LE. 0) PRINT 1235,NTOT,ITY,JTY
 1235 FORMAT('0*GHEISH* 1235 NTOT,ITY,JTY = ',3(I10,1X))
      IF (JTY .LE. 0) STOP
*     ITRT=Q(JTY+6)
      PLX=EVE(INDEX+9)
      PLY=EVE(INDEX+10)
      PLZ=EVE(INDEX+11)
      ELT=SQRT(PLX*PLX+PLY*PLY+PLZ*PLZ+Q(JTY+7)**2)
C
C --- ADD PARTICLE TO THE STACK IF STACK NOT YET FULL ---
      IF (NGKINE .GE. MXGKIN) THEN
          WRITE(CHMAIL,1236) NTOT, L
 1236     FORMAT(' *** GHEISH: ',I9,' particle produced but only ',
     +           I9,' put on the GEANT stack!')
          CALL GMAIL(1,1)
          GO TO 9999
      ENDIF
      NGKINE=NGKINE+1
      GKIN(1,NGKINE)=PLX
      GKIN(2,NGKINE)=PLY
      GKIN(3,NGKINE)=PLZ
      GKIN(4,NGKINE)=ELT
      GKIN(5,NGKINE)=ITY
      TOFD(NGKINE)=EVE(INDEX+8)*0.5E-10
      GPOS(1,NGKINE) = VECT(1)
      GPOS(2,NGKINE) = VECT(2)
      GPOS(3,NGKINE) = VECT(3)
C
      IF (NPRT(9)) PRINT 1006,ITY,NGKINE,L,(EVE(INDEX+J),J=1,12)
 1006 FORMAT(' *GHEISH* GEANT PART. ',I3,' ALSO PUT ONTO STACK AT',
     $ ' POS. ',I3/
     $ ' EVE(',I2,') = '/12(1H ,12X,G12.5/))
C
 61   CONTINUE
C
 9999 CONTINUE
C --- LIMIT THE VALUE OF NGKINE IN CASE OF OVERFLOW ---
      NGKINE=MIN(NGKINE,MXGKIN)
      END
